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Abstract 

We study the morphological evolution of strained islands in growing crystal films by use of a 
continuum description including wetting, elasticity and deposition flux. Wetting breaks transla- 
tional invariance, allowing the flux to tune different nonlinear regimes. Increasing the flux, we find 
first an annealing-like dynamics, then a slower but non-conventional ripening followed by a steady 
regime, while the island density continuously increases. The islands develop spatial correlations 
and ordering with a narrow two-peaked distance distribution and ridge-like clusters of islands for 
high flux. 

PACS numbers: 68.55.-a, 81.15.Aa, 68.35. Ct 
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Non-equilibrium crystal growth is questioning many fundamental and experimental issues 



in particular in the domain of self-organization of nanostructures 



mm 



5j . For example, 



quantum dots has lead to numerous applications Q| such as photovoltaic cells, memory 
storage or light emission. However, the different scenario governing island formation are 
still challenged as regards their density, size distribution and spatial ordering. We focus 
here on the properties of islands produced in heteroepitaxy. 

When a film is coherently deposited on a substrate with a lattice mismatch, it experiences 
an elastic stress that can be relieved by a morphological evolution. For strong enough 
mismatch, islands are nucleated in an abrupt two to three dimensions transition. However, 
for intermediate mismatch, the evolution begins by surface diffusion with a nucleationless 
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Asaro-Tiller-Grinfel'd instability [111. |12||. Contrarily to the evolution in thick films pjj], no 
dislocations are generated in thin films where instead the ripples transform smoothly into 



islands separated by a wetting layer {qI. This scenario requires ful 
nonlinear regime which involves nonlocal elastic interactions jisl . 
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understanding of the 



16 



17l | . Crucial for 



possible application, the resulting spatial order depends strongly on the growth dynamics. 
Hence, the question we address in this Letter is the influence of the dynamics on the spatial 
organization of the islands. To tackle this problem, we use a continuum modelization of 
the crystal accounting both for wetting and elastic interactions. Wetting, which breaks the 
translational invariance in the growth direction, introduces a significant flux dependence. 
Consequently, we find different nonlinear regimes as the flux F increases. Departing from 
the near annealing case, island ordering arises first as clusters form. Then, for high enough 
fluxes, the island density is frozen and the dot density and ordering are maximum. 

The dynamics of a surface during crystal growth involves different mechanisms such as 
diffusion, attachment or relaxation. The evolution of its interface z = H{r,t) at position 
r = (x, y) with time t can be written as 

flFf 

— = V[H] + ajF, (1) 

where V[H] is dictated by the predominant mechanisms at stake while aj is the film lattice 
parameter. In homoepitaxy, symmetry constraints enforce V to depend only on the slope 
of H [5] so that F disappears in the Galilean transformation H{r,t) = Ft + h{r,t). This 
invariance is violated in heteroepitaxy when a film is coherently deposited on a substrate. 
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We consider a film evolving in the Stransky-Krastanov mode, with a dynamics due to surface 
diffusion induced by chemical potential gradients. Hence, 



V[H] = D^l + \VH\'^V's^^, (2) 

where D is a diffusion coefficient, V5, the surface gradient, and /x, the surface chemical poten- 
tial given by the functional derivative /i = 5(JF^^+JF'^) / where = J dr'y{H) a/1 + |Vi?P 
and JF^' is the elastic free energy. Wetting is embedded in the if-dependence of the surface 
energy 7, which describes the chan ge w ith H of the local environment of a particle when the 
film/substrate interface is present [18'], and which precisely breaks translational invariance 
in the z direction. Note that neither alloying nor anisotropy are considered here, which can 



prove significant in some systems 19|, |20 1 . 



To fix scales, we choose to depict a Geo.75Sio.25-like film deposited on a Si substrate with 
a reference lattice with a substrate (film) lattice parameter = 0.27 nm (aj = l.Olas). 
Surface diffusion is given by -D = exp[— _E'(i/A;sT]aVA;sT with Ed = 0.83 eV and Dg = 
8.45 10-,„Vs. see e.g. a„d the working temperature is TOO^C. The film surface euergy 
is 7/ = 1.3 J/m^ and we extrapolate ab-initio calculations for Si/Ge systems [21] by consid- 
ering 7(/i)=7j[l + Cwexp{—h/6)] where 5 = aj and c^ = 0.09. The characteristic length and 
time scales are then Iq = £^ /2{1 + ^'/)7/ and 10 = 1^/ D'jf, where £^ = Ef {aj — as)"^ /a^^l — Uf) 
is an elastic energy density involving the film Poisson ratio and Young modulus Uf and Ef, 
whose values are /o = 27nm and to = 25 s. 

Mechanical equilibrium is supposed to be achieved on time scales faster than the system 
evolution, which enforce the Lame equations dgapg = to be valid in the film and substrate. 
The stress tensor (Jpq IS Sb linear function of the strain tensor Cpg = ^{dqUp + dpUq) — Cp^ where 
u is the displacement with respect to the substrate reference state. The reference strain is 
epg = {l — af / as) 5pq{5px + 5qy) in the film and otherwise, where 5p^q is the Kronecker symbol 
with p,q = x,y,z. To simplify calculations, we consider isotropic elasticity and identical 
film and substrate elastic constants as coherent epitaxy involves similar materials. This 
approximation neglects higher order terms and ensures the translational invariance with z 
of the elastic energy; consequently, in the Galilean frame h = H — h{t), the dependence on 
h{t) will appear exclusively through the wetting term. The Lame equations are then solved 
thanks to Fourier transforms with respect to r [id] with the following boundary conditions : 
the film/ substrate interface is coherent with continuous stresses, and the film/vacuum surface 
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is free, with a negligible surface stress. To solve the latter condition, we use here the small- 
slope approximation amenable for arbitrary deposited thicknesses, contrarily to the thin 



film approximation used in [16|. Writing h the spatial average of /i, we suppose a shallow 
modulation where h — his small compared to the lateral characteristic length of order Iq, and 
get displacements and /i^', the elastic energy density ^(Jpgepg computed at the film surface, 
up to second order in this small parameter. Considering h = afFt, the second term in the 
r.h.s. of Eq. ([1]) cancels out, so that 

r)h r / - \ r p-{h+Ft)/S 

? = A I - f 1 + c^e-('^+^*)/^) Ah - ^ " 

-Hiiih) +2hAh+\Vhf 

{h)e,,mki{h)] , (3) 

in units of /q and to and with F = afFto/lQ. The long-range elastic interactions enforce the 
non-analytic behavior of the operator 

n^Ah]=:F-'{{hk,/\k\)j^[h]}, (4) 

defined with Fourier transforms JF over r and wavevector k, where the indices run over 
x,y. The nonlinear terms are given with Oijij = 1, Oajj = —Oijji = Vf when i ^ j and 
otherwise. It is worthwhile to realize that a lone term hAh in the r.h.s. of (131) would be 
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ruled out by symmetry considerations [5] but is allowed here when put in balance with the 
nonlocal nonlinear terms. Indeed, taking into account both local and nonlocal terms ensure 
the invariance of the elastic energy under the transformation h—>-h + h{t) for arbitrary h{t). 

We performed numerical simulations of Eq. using a pseudo-spectral method. We 
consider an initial roughness of amplitude 1 monolayer (ML) given by a random profile. 
The initial film height is ho{t = 0) = 7 ML just below the elastic instability threshold he — 
8.2 ML defined below. Similarly to the annealing case liol], the simulations reveal that 
the combination of wetting and nonlinear nonlocal elastic terms prevents the finite-time 
singularities, which in thick films, lead to dislocations In addition, we find that the 

system evolution depends strongly on the deposition flux. Different curves of the roughness 
w = [{h^) — {h)'^Y^'^ as function of the deposited height h^ = Ft are depicted in Fig. m for 
different fluxes. The roughness first increases exponentially in a linear-like dynamics. After a 
first inflexion point, the system enters a first nonlinear stage where w^t, and after a second 




FIG. 1: Simulation results of the roughness evolution with the deposited height for, from left to 
right, F = 10-^510-^,10-^510-^10-^210-^ 510-2,10-iML/s. 

one, displays a faster than linear roughness increase. The deposited heights corresponding 
to these inflexions increase linearly with the flux, though the different curves cannot be 
rescaled on a single one, signaling different ripening mechanisms depending on the flux. 

To evaluate the time of emergence of the islands, we define a dynamical critical height 
hs{F) after which w is greater than Wg = 3 ML. It is well fitted by an affine law hg = 
10+1200-F (ML) with a limit at low flux greater than he due to the threshold present in the 
definition of h^. The values of /i, obtained here for T = 700° C differ from the smaller values 



of the apparent critical thickness derived with the analysis of 121]. The latter is defined via 
a comparison of relative growth rates, whereas hg is defined here directly from the roughness 
which is independent of the reference frame. The dynamical critical height hg is related to 
the growth in the linear regime and can be roughly estimated. In the linear approximation, 
h small, the evolution of ([3]) can be reformulated in Fourier space as 

h{k; t) = h{k] to) exp ds a (^k; Jiq + Fs)] , (5) 

where 

a(k; h) = -fc'c^e-^/V'^' + I'^l' - k\l + c^e"'^/'^). (6) 

By definition of he, o is negative for h < he and otherwise positive in a /i-dependent k- 
interval. To estimate hg, we consider an initial undulation with wavevector k^, and an initial 
amplitude hi^ = 1 ML. We can solve exactly for hl{r; t) and the critical deposited height and 
find hl{F) = ^ + 6W{-3c^{l + k'^6^) exp{-^/5)/kl6^), with the product-log function W 
and the length ^ = ho + 4c^(l + kl6'^)e~-^°/^ /kj + F\og{wg/hi^o)/kt- Considering A;* = 0.35, 
the solution is well fitted by hl{F) = 10+1300F (ML) which is a rough estimate of hg for the 



FIG. 2: Typical island configurations in the nonlinear regime for F = (a) 10 (b) 10 ^, (c) 
10~^ML/s in a 128/oxl28/o system. 

full equation at low flux. Similarly to the numerical results, h* increases nearly linearly with 
F in the regime of parameters studied here. Also, the limit at low fluxes is greater than the 
instability critical height he due to the threshold in the definition of hg. These results could 
be confronted to experiments investigating the appearance of islands as function of the ratio 
D/F. 

Above the critical height, the system enters the nonlinear regime characterized by islands 
surrounded by a wetting layer growing both by deposition and coarsening, see Fig. [2j The 
island density p is then plotted in Fig. [3] for different fluxes. The islands appear all the 
later that the flux is high, as a result of the competition between the instability time scale 
and the deposition growth. In addition, the dynamics depends strongly on the flux. At low 
flux, the island density is convex, similarly to the annealing case, when the system has time 
to coarsen by surface diffusion. However, in the intermediate regime for F in between Fi 
and F2 with Fi~10~^ML/s and F2 ~5 10~^ ML/s, the density evolution becomes concave, 
whereas in the steady regime, F>F2, the island density is constant in a large time interval. 
The values of Fi and F2 depend on the temperature and the details of the wetting that are 
set here with typical coefficients. In all regimes, the island density, together with their areal 
coverage, at a given deposited height, is an increasing function of the flux and saturates at 
a value limited by the linear instability wavelength. This observation indicates a route for 
controlling the island density. We find also that, in the intermediate and steady regimes, 
the areal coverage first increases with the deposited height even in the nonlinear stage. 

The change of dynamics from a convex island density evolution in the annealing-like case 
to a concave one for growing films is similar to the results of the experiments by Floro and co- 
workers |lo| where the negativeness of (Pp/dt^ is signaling a non standard ripening. Indeed, 
a typical coarsening dynamics tends to decrease its driving force resulting in a damped 
evolution. A faster dynamics can nevertheless be apprehended within a mean-field approach 
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FIG. 3: Island density as function of the deposited height with, from left to right, the fluxes of 

Fig.m 

describing the island size distribution evolution due to chemical potentials accounting for 

n 

elastic interactions [10]. These long range interactions are precisely the central ingredient 
of ([3]) which can lead to such a dynamics when coupled with deposition growth and our 
simulations exhibit a regime where this concave evolution is to be expected. Moreover, in 
all regimes, the dot density evolution is slowed down even when it is given as a function 
of the instability time scale instead of the deposited height. Hence, even if a wetting layer 
has developed and allows surface diffusion, the stabilization of the islands results from an 
effective weakening of surface diffusion currents when deposition is present. 

The different nonlinear regimes are moreover characterized by different spatial organiza- 
tion, see Fig. [2l At low flux, the islands do not exhibit strong spatial correlations, similarly 
to the annealing case [l^, while for higher flux, island decimation, driven by surface dif- 
fusion, leads to the emergence of island clusters. These clusters, absent in the annealing 
case, involve more and more islands as F increases (Fig. [2]). To quantify these correlations, 
we first assign an area mass center to each island and construct a Delaunay triangulation 
for this set of points. The typical nearest neighbor distance histograms in the nonlinear 
regimes are then plotted in Fig. HI At low flux, the histogram displays a broad bell shape 
with a modulation reminiscent of the instability initial stages, where, beside the first peak 
near the initial wavelength A ^ 201q, the other peaks arise as coarsening is fully developed. 
For intermediate flux, the histogram is significantly narrowed with a main peak near A, and 
other peaks which positions increase slightly with time as island decimation occurs. In the 
steady high flux regime. Fig. Ht), ordering is maximal and is described by a narrow dis- 
tance distribution with two peaks related to the typical distance between two islands in a 
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FIG. 4: Island distance distribution corresponding to Fig. [2] 
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FIG. 5: The kinetic phase diagram as function of the deposited height h(i = Ft and the flux. For 
a given flux, growth proceeds mainly layer by layer up to the dynamical critical height hs (solid 
line). Above, for F<Fi, surface diffusion is efficient and an annealing-like ripening proceeds. For 
Fi < F < F2, a non-conventional coarsening is at stake where the island density decreases faster 
than linearly. Finally, for F > F2, ripening is frozen by the deposition growth, the island density 
stays constant and highly correlated islands tend to form ridge- like structures. 



cluster and between clusters. The first peak arises at a position significantly lower than the 
instability wavelength A, signaling a narrowing of the distance between mass centers during 
clustering at high flux, while the second peak sticks at A. In this case, a particular ordering 
is observed in Fig. [Jt) as islands self-organize in ridge-like patterns which result from elastic 
interactions and stabilization of the nonlinear regime by the flux. Finally, note that the 
histograms of the Voronoi tessellation do not give relevant information here, contrarily to 
systems with island nucleation. These results are summarized by the kinetic phase diagram 
in Fig. [51 

As a conclusion, we studied the influence of a deposition flux on the growth dynamics 
leading to island formation. We used a continuum description of the elastic instability 
accounting for wetting and nonlocal elasticity effects where the flux arises in a nontrivial 
way. We found three different nonlinear regimes, depending on the flux and characterized 
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by different spatial ordering and dynamics. At low flux, an annealing-like dynamics is at 
stake. Increasing the flux, we find spatial correlations where dots gather in clusters, together 
with a non-conventional ripening. For high flux, ripening is frozen as surface diffusion effects 
are inhibited by deposition growth and the dot density is maximum. The nearest neighbor 
distance distribution exhibits a first peak linked to the distance between two islands in 
a cluster which decreases, and a second one related to the cluster distance. These results 
indicate a way to tune different spatial ordering and are calling for experimental examination. 
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